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Abstract 

We examine a model system where attractors may consist of a heteroclinic cycle 
between chaotic sets; this 'cycling chaos' manifests itself as trajectories that spend 
increasingly long periods lingering near chaotic invariant sets interspersed with short 
transitions between neighbourhoods of these sets. This behaviour can be robust (i.e., 
structurally stable) for systems with symmetries and provides robust examples of 
non-ergodic attractors in such systems; we examine bifurcations of this state. 

We discuss a scenario where an attracting cycling chaotic state is created at a 
blowout bifurcation of a chaotic attractor in an invariant subspace. This is a novel 
scenario for the blowout bifurcation and causes us to introduce the idea of set 
supercriticality to recognise such bifurcations. The robust cycling chaotic state can 
be followed to a point where it loses stability at a resonance bifurcation and creates 
a series of large period attractors. 

The model we consider is a 9th order truncated ordinary differential equation 
(ODE) model of three-dimensional incompressible convection in a plane layer of 
conducting fluid subjected to a vertical magnetic field and a vertical temperature 
gradient. Symmetries of the model lead to the existence of invariant subspaces for 
the dynamics; in particular there are invariant subspaces that correspond to regimes 
of two-dimensional flows. Stable two-dimensional chaotic flow can go unstable to 
three-dimensional flow via the cross-roll instability. We show how the bifurcations 
mentioned above can be located by examination of various transverse Liapunov 
exponents. We also consider a reduction of the ODE to a map and demonstrate 
that the same behaviour can be found in the corresponding map. This allows us to 
describe and predict a number of observed transitions in these models. 
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1 Introduction 

There has been a lot of recent interest in the chaotic dynamics of nonlinear 
systems that possess invariant subspaces; a number of quite subtle dynamical 
effects come to light in examining the interaction of attractors with the in- 
variant subspaces. Moreover a good understanding of these effects is essential 
in interpreting and predicting dynamics of simulations and experiments where 
the presence of discrete spatial symmetries implies the existence of invariant 
subspaces. 

A fundamental bifurcation in such a setting is the blowout bifurcation [16] 
where a chaotic attractor within an invariant subspace loses stability to per- 
turbations transverse to the invariant subspace. In doing so, can either create 
a nearby 'branch' of chaotic attractors or lose stability altogether. 

In contrast with this, the presence of invariant subspaces can lead to the exis- 
tence of what have been called robust heteroclinic cycles [11] between equilib- 
ria, that is, heteroclinic cycles that are persistent under small perturbations 
that preserve the symmetry. These cycles may or may not be attracting [13]. 
Recently it has been recognised that cycles to more complicated invariant sets 
can also occur robustly in symmetric systems, in particular to chaotic invariant 
sets; this behaviour was named 'cycling chaos' in a recent paper by Dellnitz 
et al. [9] and has been further investigated by Field [10] and Ashwin [2]. 

In this paper we find there is a connection between these dynamical properties; 
we show a scenario where a blowout bifurcation creates an attracting 'cycling 
chaotic' state in a bifurcation that is analogous to a saddle-node homoclinic 
bifurcation with equilibria replaced by chaotic invariant sets. We also inves- 
tigate how the attracting cycling chaotic state that is created in the blowout 
bifurcation loses stability at a resonance of Liapunov exponents. (A resonance 
bifurcation in its simplest form occurs when a homoclinic cycle to an equi- 
librium loses attractiveness; this occurs when the real parts of eigenvalues of 
the linearisation become equal in magnitude [7].) In spite of the system being 
neither a skew product nor being a homoclinic cycle to a chaotic set as in [2] 
we see similar behaviour and can predict the loss of stability by looking at a 
rational combination of Liapunov exponents. 

We find this scenario of a blowout bifurcation to cycling chaos is a mech- 
anism for transition from stable two-dimensional to fully three-dimensional 
magnetoconvection. The model we study is a Galerkin truncation for magneto- 
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convection in a region with square plan, subject to periodic boundary condi- 
tions on vertical boundaries and idealised boundary conditions on horizontal 
boundaries. Phenomenologically speaking we see a change from a chaotically 
varying two-dimensional flow (with trivial dependence on the third coordinate, 
and which comes arbitrarily close to a trivial conduction state) to an attracting 
state where trajectories spend increasingly long times near one of two sym- 
metrically related two-dimensional flows interspersed with short transients. 
We explain and investigate this transition in terms of a blowout bifurcation 
of a chaotic attractor in an invariant subspace. 

In the paper of Ott and Sommerer [16] that coined the phrase 'blowout bi- 
furcation', two scenarios are identified. Either the blowout was supercritical 
in which case it leads to an on- off intermittent state [17], or it is subcritical 
and there is no nearby attractor after the bifurcation. We find an additional 
robust possibility for bifurcation at blowout. 

Near this transition the three-dimensional flow patterns show characteristics 
of intermittent cycling between two symmetrically related 'laminar' states 
corresponding to two-dimensional flows, but the time spent near the laminar 
state is, on average, infinite. This suggests that the blowout is supercritical, 
but in a weaker sense which we make precise. Namely, we say a blowout is set 
supercritical if there is a branch of chaotic attractors after the blowout whose 
limit contains the attractor in the invariant subspaces before the blowout. In 
particular there may be other invariant sets contained in this limit and so any 
natural measures on the bifurcating branch of attractors (if they exist) need 
not limit to the natural measure of the system on the invariant subspace. 

We also show that the attractors corresponding to two-dimensional flows are 
not Liapunov stable, but are Milnor attractors near the transition to three 
dimensions, and so in particular we expect the presence of noise to destabilise 
two-dimensional attractors near blowout by a bubbling type of mechanism [4]. 

We find in our example that the state of cycling chaos is attracting once it 
has been created: trajectories cycle between neighbourhoods of the chaotic 
sets within the invariant subspaces, and the time between switches from one 
neighbourhood to the next increases geometrically as trajectories get closer 
and closer to the invariant subspaces. By estimating the rate of increase of 
switching times, we are able to show that cycling chaos ceases to be attracting 
in a resonance bifurcation. One remarkable aspect of this study is that we are 
able to predict the parameter values at which the blowout bifurcation and the 
resonance occur, requiring only a single numerical average over the chaotic set 
within the invariant subspace. 

The paper is organised as follows: in section 2 we introduce the ODE model for 
magnetoconvection, discuss its symmetries and corresponding invariant sub- 
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spaces. This is followed by a description of the creation, persistence and loss 
of stability of the cycling chaos on varying a parameter in numerical simula- 
tions in Section 3. Section 4 shows how one can, under certain assumptions, 
derive a map model of the dynamics of the ODE that has the same dynamical 
behaviour. Section 5 is a theoretical analysis of the blowout bifurcation that 
creates the cycling chaotic attractor and is followed in Section 6 by a theo- 
retical analysis of its loss of stability. Finally, Section 7 discusses some of the 
implications of this work on the chaotic dynamics of symmetric systems. 



2 An ODE model for magnetoconvection 

The model we study is an ODE on R 9 described by the following equations 



x 


= /ix + x 9 - x 2 x 1 - PvqX , 




= —VX\ + X§X2i 




= —ax2 — crQa + 7x0^1, 


yo 


= Wo + yoO - 2/22/1 - Pxly , 


yi 


= -vyi + y y 2 , 


h 


= -ay 2 - oQb + 72/oJ/i, 


a 


= ((x 2 - a), 


b 


= C(lfa-&), 


6 


= -9-x 2 -yl 



These ODEs have been derived as an asymptotic limit of a model of three- 
dimensional incompressible convection in a plane layer, with an imposed verti- 
cal magnetic field; for further details and details of its derivation, see [14,19,20]. 
In the context of this model, x and y represent the amplitudes of convective 
rolls with their axes aligned in the y and x (horizontal) directions respectively, 
x\ and y\ represent modes that cause the rolls to tilt, and x 2 and y 2 represent 
shear across the layer in the x and y directions. The modes a and b represent 
the horizontal component of the magnetic field in the x and y directions, and 
6 represents the horizontally averaged temperature. 

The model has are five primary parameters: /i is proportional to the imposed 
temperature difference across the layer, with /i — at the initial bifurcation to 
convection; (3 is related to the horizontal spatial periodicity length, but is an 
arbitrary small parameter in the model of [14,19]; a and ( are dimensionless 
viscous and magnetic diffusion coefficients; and Q is proportional to the square 
of the imposed magnetic field. Note that a, ( and Q are scaled by factors of 
4, 4 and 7r 2 from their usage in [14,19,20]. Two secondary parameters that we 
use are 7 = 3(1 + 4cr)/16er and v = (9cr/(l + 4cr)) — ji. In the parameter regime 
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Table 1 

Selected fixed point subspaces S of the action of G on R 9 together with name, 
representative point and dimension of S. There are others (e.g. (0, x, 0, 0, 0, 0, 0, 0, *)) 
but these are not important for the dynamics we discuss here. 



s 


Name 


Representative point 


dimS" 


F 


Full symmetry 


(0,0,0,0,0,0,0,0,*) 


1 


Rx 


x-rolls 


(x, 0,0, 0,0,0,0,0,*) 


2 


Ry 


y-rolls 


(0,0,0,3!, 0,0,0,0,*) 


2 


D+ 


+ diagonal 


(£,0,0,2,0,0,0,0,*) 


2 


D 




fr -r o o n f) 


2 


Rxy 


Mixed modes 


(ar, 0,0, j/, 0,0,0,0,*) 


3 


P X 


x-rolls + shear 


(x, y, z, 0, 0, 0, a, 0, *) 


5 


Py 


y-rolls + shear 


(0, 0, 0, z, y, 0, a, *) 


5 


Qx 


x-rolls + shear + crossrolls 


(x, y, z, w, 0, 0, a, 0, *) 


6 


Qy 


y-rolls + shear + crossrolls 


(w, 0, 0, x, y, z, 0, a, *) 


6 


T 


No symmetry 


(u, v, w, x, y, z, a, b, t) 


9 



of interest, all parameters are non-negative. 



2. 1 Symmetries of the model 

Consider D 4 -i-T 2 acting on the plane with unit cell [0, 2n) 2 in the usual way, 
with the torus T 2 acting by translations on the plane, and D 4 by reflections 
in the axes and rotation through ir/2. We define the following group elements 



k x : 


reflection through x = 





(x,y) h 




x,y), 


K 'x : 


reflection through x = 


tt/2 


(x,y) h 


-> (ir 


-x,y), 


K y : 


reflection through y = 





(x,y) h 


-> (x 


-y), 


n' y . 


reflection through y = 


tt/2 


(x,y) h 


■+ (x 


n-y), 


P- 


rotation about (x, y) = 


= (7T,7r) 


(x,y) h 


-> (2ir-y,x), 




translation 




(x,y) h 


-> (x 





(2) 



Note that p, t^, v ) and any reflection k can be used to generate the group 
D 4 +T 2 . 
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Fig. 1. A portion of the isotropy lattice for the action of G on R 9 under which (1) 
is equi variant. We have shown fixed point subspaces of some conjugate subgroups 
separately for clarity. The isotropies of P x and P y are the smallest isotropics that 
physically involve only two-dimensional effects. 

We consider the subgroup 

G = (k x , k' x , K y , k' ,p). (3) 



Since G contains the subgroup (Z 2 ) 2 , generated by k x k' x and K y K y , of trans- 
lations T 2 , it follows that G is isomorphic to a semidirect product D 4 +(Z 2 ) 2 
(\G\ = 32). The ODE (1) is equivariant under the group G of symmetries 
acting on R 9 by 

k x (x , x 1 ,x 2 , y , Vi, 2/2, a, b, 9) = (x , -xi, -x 2 , y , yi, y 2 , -a, 6, &), 
k' x (x , x 1: x 2 , y , yi, y 2 , a, b, 6) = (-x , x 1: -x 2 , y , y u y 2 , -a, b, 9), (4) 
p(x Q , x u x 2 , y , yi, y 2 , a, b, 6) = (y , -y u -y 2 , x , x u x 2 , -b, a, 6). 

This action gives rise to a number of isotropy types, shown in Table 1. Figure 1 
gives a partial isotropy lattice for this group action. Dynamics in F always 
decays to the trivial equilibrium point, corresponding to the absence of con- 
vection. We refer to dynamics in P x and P y as two-dimensional, since these 
correspond to two-dimensional convection in the original problem (though 
dimPj, is 5). Dynamics in R x and R y corresponds to mirror symmetric two- 
dimensional rolls with their axes aligned along the y and x directions; we 
refer to equilibrium points in these subspaces as x-rolls and y-rolls respec- 
tively. In P x and P y , convection is two-dimensional but not mirror symmetric, 
and is referred to as tilted rolls. Dynamics in Q x and Q y corresponds to three- 
dimensional convection that is still invariant under one mirror reflection (tilted 
rolls with a cross- roll component). Otherwise we say the dynamics is fully three 
dimensional. 
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In a slight break from convention we say the fixed point subspaces S as having 
isotropy subgroup Iso(jS') rather than considering the isotropy subgroups as 
the fundamental objects. 



3 Numerical simulations of the ODEs 

We present numerical simulations of the ODEs that demonstrate two aspects 
of cycling chaos that we seek to explain: how cycling chaos can be created in 
a blowout bifurcation, and how cycling chaos can cease to be attracting. 

We concentrate on parameter values that are known to have Lorenz-like chaotic 
dynamics within P x and P y : 

^ = 0.1655, Q = l/n 2 , a = 0.125 and C = 0.05 (5) 

(and hence v = 0.5845 and 7 = 2.25), although we note that qualitatively 
similar attractors are found for a large proportion of nearby parameter values. 
These parameter values correspond to those in Figures 15(c) and 20(a) in 
[20]. The numerical method used was a Bulirsch-Stoer adaptive integrator 
[18], with a tolerance for the relative error set to 10~ 12 for each step. 

We use the parameter (3 as a normal parameter (see [5]) for the dynamics 
in P x and P y \ that is, it controls instabilities transverse to P x and P y in the 
directions Q x and Q y without altering the dynamics inside P x or P y . 

3.1 Cycling chaos 

In Figure 2 we show a typical example of timeseries when there is attracting 
cycling chaos (with (3 = 2.0). The system starts with xq oscillating chaotically 
while yo is quiescent, switches to a state where yo oscillates chaotically while xq 
is quiescent, and so on. A more careful examination reveals that after a switch 
the trajectory remains close to a fixed point in R x or R y for an increasing length 
of time. Physically, this corresponds to chaotic two-dimensional convection 
that switches between rolls aligned in the x and y directions. Figure 2(c) shows 
the chaotic trajectories projected onto the (2:0, ^2) plane, while (d) illustrates 
switching between P x (the 'horizontal' plane) and P y (the 'vertical' plane). 

Note how the chaotic behaviour in Figure 2(a) and (b) repeats: trajectories 
spend longer and longer near the unstable manifolds of the rc-roll and y-roll 
equilibrium points and take longer and longer between each switch. This is 
illustrated further in Figure 3, which shows intersections of a trajectory with 
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Fig. 2. Numerical solutions of the model ODEs with [i = 0.1655, Q = l/n 2 , 
a = 0.125, £ = 0.05 and (3 = 2.0. (a) xo against time; (b) j/o against time; (c) xo 
against x 2 ; (d) Xq, and |x 2 | + I2/2I i n perspective. The crosses in (c) and (d) rep- 
resent the x-roll, y-roll and trivial equilibrium points in R x , R y and F. Note that 
the chaotic attractor is close to a heteroclinic cycle that connects the x-rolls, the 
y-rolls and the trivial fixed point. 
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Fig. 3. Numerical solutions of the model ODEs: cycling chaos with (a) (3 = 2.00, 
(b) (3 = 1.80, (c) (3 = 1.65 and (d) (3 = 1.63, showing the values of |xo| at which the 
trajectory intersects the surface defined by \x 2 \ + \y2\ = 0.01. The rising exponential 
growth phases correspond to the system switching from chaos in P y to chaos in 
P x . The time between switches increases as the system approaches the attracting 
cycling chaos, but the rate of increase depends on (3. 

the Poincare surface \x 2 \ + 1 2/2 1 = 0.01 close to the trivial solution. There 
are two phases evident in the cycle: the order one chaotic behaviour of xo 
near P x (while i/q grows exponentially), and the exponential growth of |xo| 
as the trajectory moves away from P y (while y behaves chaotically). The 
time between switches increases monotonically and the rate of increase varies 
with (3, the normal parameter. In these numerical simulations, the switching 
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Fig. 4. A schematic illustration of the location of the robust cycle relative to the 
invariant subspaces forced by symmetry. The cycle is between the chaotic invariant 
sets A x and A y (within P x and P y ) and two fixed points contained in R x and R y . 
The cycle is robust to G-equivariant perturbations that fix the dynamics in P x and 
P y . The intersection of P x and P y at the trivial solution is not shown, although A x 
and A y do in fact intersect there. 

time saturates when certain components of the solution come close to the 
machine accuracy of the computations (about 1CT 323 for double precision). 

We argue that this is evidence for attracting cycling chaos: trajectories ap- 
proach a structurally stable heteroclinic cycle between chaotic sets. In Fig- 
ure 4 we show a schematic picture of the heteroclinic cycle. We recall that the 
fixed point spaces R x and R y are 2 dimensional, R xy is 3 dimensional, P x and 
P y are 5 dimensional and Q x and Q y are 6 dimensional invariant subspaces in 
the 9 dimensional phase space. The system starts near the x-roll equilibrium 
point in R x , which is unstable to shear (#2); the parameters are such that the 
unstable manifold of rr-rolls comes close to the trivial solution and returns to a 
neighbourhood of :r-rolls. This global near-connection within P x is the source 
of the chaotic behaviour. We refer to the chaotic sets in P x and P y as A x and 
A y respectively; these contain the relevant roll equilibrium points and they 
both contain the trivial solution, so there are structurally stable connections 
from the trivial solution to the roll equilibrium points and from those to A x 
and Ay. Within Q x , A x is unstable to cross-rolls (y ) since the trivial solu- 
tion is equally unstable in the x and y directions. Eventually y grows large 
enough that there is a switch to the y-roll equilibrium point in R y , at which 
point i/2 starts to grow. Thus the cycle connects invariant sets in the following 
fixed point subspaces: 

. . . > R x > P x > Q x ^ Ry ^ P y ^ Qy ^ • • • (6) 

between the equilibrium points in R x and R y , and A x and A y within P x and P y , 
with the structurally stable connections needed to complete the cycle from A x 
and A y to the ?/-roll and a;-roll equilibrium points lying within Q x and Q y . 
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Note that our scenario is certainly a simplification of the full set of hetero- 
clinic connections; there are other fixed points contained in the closure of the 
smallest attracting invariant set, in particular the origin is contained within 
the cycle. 



3.2 Blowout 



We turn now to the question of how the cycling chaos is created. With (3 = 3.5 
we find attracting two-dimensional chaos (Figure 5a, b), which loses stability 
around (3 = 3.47 (c) in a blowout bifurcation, and for (3 = 3.40 (d), there is 
exponential growth away from P x into Q x . Within Q x , the y-roll equilibrium 
points are sinks, establishing the structurally stable connection from the chaos 
in P x to the equilibrium point in R y , and hence the creation of cycling chaos. 



3.3 Resonance 

As illustrated in Figures 3, as (3 decreases towards 1.63, trajectories spend 
progressively longer in each visit to P x and P y before switching to the conju- 
gate chaotic invariant set. Eventually, trajectories come arbitrarily close to the 
invariant subspaces P x and P y (limited only by machine accuracy in the nu- 
merical simulations). Figure 6 shows how the time intervals between switches 
between A x and A y increases as the system approaches this heteroclinic cycle, 
and how the rate of approach to the heteroclinic cycle decreases as (3 ap- 
proaches 1.62. The intervals between switches grow by a factor of about 1.4, 

1.2 and 1.1 per switch for (3 = 2.00, (3 = 1.80 and (3 = 1.65 respectively, and 

1.03 and 1.02 for (3 = 1.64 and 1.63. By (3 ~ 1.62, the heteroclinic cycle is 
no longer attracting, and for (3 = 1.61 and 1.60, the system settles down to 
periodic behaviour that is bounded away from P x and P y , though the periodic 
orbits are actually quite close to these invariant subspaces. For these calcula- 
tions, we imposed a cut-off of 10 -100 : the calculation ceased once any variable 
became smaller than this. 

In the next section, we derive a map that allows us to compute longer tra- 
jectories more accurately, and we demonstrate cycling chaos, its creation in a 
blowout bifurcation and its loss of attractiveness using this map. We analyse 
the blowout bifurcation in section 5, and argue in section 6 that the cycling 
chaos created in that bifurcation ceases to be attractive at a resonance. 
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Fig. 5. Numerical solutions of the model ODEs. The top panel (a) shows the chaotic 
time dependence of xq, independent of (3, and the lower three panels show the 
blowout bifurcation: (b) (3 = 3.50, before the blowout bifurcation (solutions de- 
cay to two-dimensional chaos); (c) (3 = 3.47 near the blowout bifurcation; and 
(d) (3 = 3.40 after the blowout bifurcation (solutions grow to cycling chaos). The 
straight lines show the average growth or decay rates predicted by the Liapunov 
exponents calculated using (27). 
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Fig. 6. Time intervals between switches between P x and P y with /? in the range 
1.60-1.64. The intervals between switches increase by a factor of about 1.03 and 1.02 
per switch for (3 = 1.64 and 1.63 respectively. The switching times stop increasing 
once the values of variables become so small that the cubic terms in the ODEs 
cannot be represented accurately with double-precision numbers (and calculations 
were terminated once any variable went below 10 -100 ). For (3 = 1.61 and 1.60, the 
system approaches a periodic orbit and the intervals between switches goes to a 
constant (though in the case of (3 = 1.61, the periodic orbit is only achieved after 
1500 switches). 

4 Reduction to a map 

In this section, we discuss a map that models the behaviour of the ODEs in 
the parameter regime described above, and that helps clarify what happens 
near blowout bifurcation and the resonance of the cycling chaos. 

4-1 Derivation 

We rely on a map derived for the two-dimensional dynamics by breaking up 
the flow into pieces near and between equilibrium points [20]. Within the P x 
subspace, the leading stable eigendirection of the origin (that is, the eigendi- 
rection with negative eigenvalue closest to zero) is in the (x2, a) plane, and the 
unstable direction of the origin is along the x axis. All trajectories leaving 
the origin in that direction follow the structurally stable connection within R x 
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to one or other of the x-roll equilibrium points. The one-dimensional unsta- 
ble manifold of a;- rolls lies within P x , and the chaotic attractor is associated 
with a global bifurcation at which that unstable manifold collides with the 
origin. Near this global bifurcation, the dynamics within P x is modelled by an 
augmented Lorenz map: 

(x ,x 2 = ±1) h-> (sgn(x )(-K + d\x \ 51 ), -rr 2 ) (7) 

defined as a map from the surface of section \x 2 \ = h to itself. Details of the 
derivation are given in [20], but briefly, k is a parameter related to /i in the 
ODEs, with k = at the global bifurcation; h is a small positive constant that 
we scale to one; sgn(x) = +1 if x > and — 1 if x < 0; 5i depends on the 
ratio of various stable and unstable eigenvalues at the origin and the x-roll 
equilibrium points (with < 5i < 1); and C\ is a (negative) constant. 

It is a straightforward matter to include the effect of a small perturbation in 
the yo and y 2 directions. Near the P x subspace, y an d y 2 will grow linearly at 
a rate that depends on x . This means that we get a mapping of the form 

(x ,x 2 = ±l,J/o,J/2) i-> 

(sgn(x )(-^ + Cilxol- 51 ), -x 2 , C^okol'' 2 , C 3 2/ 2 ko|' 53 ) , (8) 

where C 2 and C3 are constants and the exponents S 2 and £3 again depend on 
ratios of eigenvalues, with 5 2 < and S 3 > 0. The exponent <5 2 is negative 
since a small value of a; implies that the trajectory spends longer near the 
origin, so yo has a longer time to grow. Similarly, near P y we get the mapping 

(x ,x 2 ,y ,y 2 = ±1) i-> 

(C 2 x |y | 52 ,C 3a : 2 | 2 / | 53 ,sgn(y )(-/ t + C 1 |yo| il ),-2/2), (9) 
as long as Xo is much smaller than y . 

These maps are valid provided that the trajectory remains close to the P x 
or P y subspaces. We model the switch from behaviour near P x to near P y 
by assuming that (8) holds provided that |rr | > |yo|j otherwise, the trajec- 
tory leaves the neighbourhood of the origin along the |/o-axis, visits a y-io\\ 
equilibrium point and returns to the surface of section \y 2 \ = 1 near the origin: 

(x ,x 2 = ±l,yo,y 2 )^(C A x \y 2 \% \ 5 \±C 5 \y 2 \ 5; >\y \ 53 , (10) 

sgn(y )(-K + C 6 \y 2 \ S6 \y \ Sl ),sgn(y 2 )') , 
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where, as above, C4, C5 and Cq are constants and 84, £5 and 5q are ratios of 
eigenvalues. Similarly, (9) holds if \x \ < \y \; otherwise the trajectory switches 
from P y to P x : 

(x ,x 2 ,y ,y2 = ±1) •-»■ (sgn(x )(-K + C 6 \x 2 \ S(i \x \ Sl ),sgn(x 2 ), (11) 

C m \x 2 \ 5 *\x \ 52 ,±C 5 \x 2 \ S5 \x \ S3 ). 

Then in the full map, as long as |xo| > \yo\, the trajectory behaves chaotically 
under (8) (near the P x subspace), with x obeying a Lorenz map and y 
growing or decaying according to the value of xq. If yo grows sufficiently that 
1 2/o I > l^o |; there is one iterate of map (10), which makes xq small as the 
system switches from P x to P y , followed by many iterates of (9), and one 
iterate of (11) as the system switches back to P x . 

The S exponents can be determined from the eigenvalues of the trivial solution 
and of rr-rolls. The relevant eigenvalues of the origin are (in the notation of 
[20]) the growth rate ji of x and yo and the slowest decay rate A^ of (x2,a), 
determined by the eigenvalue of 

[7 7 e ] < 12 > 



closest to zero. The relevant eigenvalues of x-rolls are the growth rate A + and 
the slowest decay rate A 2 of (x 1 ,X2,a), determined by the eigenvalues of 



—v ^JJi 

c -c , 



(13) 



The other important eigenvalues at the x-roll equilibrium point are the decay 
rate of y (—/3fj) and the decay rate of y 2 (A^). Then we have 



Si 



5 4 = 



-A 



Pa 

A + ' 



(1 + S 6 ) , 5 2 



-A, 



A 



+ ' 



64-1, 



r -A 2 " 
Or = 



A+ • 



(14) 



From fitting the map to trajectories of the ODEs with /3 = 3.47 (see Figure 5b), 
we find values for the constants: 



d = -0.12450, 
C 4 = 609.65, 



k= -0.014072, 
C 2 = 0.27200, C 3 = 0.16791, 
C 5 = 0.86515, C 6 = 0.0015974. 



(15) 
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The eigenvalues that do not depend on (3 are: 

A c = -0.059697, A + = 0.29606, A^ = -0.048978. (16) 

With these eigenvalues, we have 6i = 0.42038, 5 3 = 0.43344, 5 5 = 0.20164 and 
5 6 = 0.16543, while S 2 (negative) and 5 4 depend on (5. 

4-2 Iterating the map 

We seek to reproduce in the map what we have observed in the ODEs: the 
blowout bifurcation that creates cycling chaos, and the loss of attractiveness of 
cycling chaos. The first of these (Figure 7) is straightforward: typical Lorenz- 
like chaos is shown in Figure 7(a), while the change from decay towards P x to 
growth away from P x at (3 ~ 3.87 is shown in (b-d). 

Cycling chaos is found after the blowout bifurcation (Figure 8a): the system 
switches between chaos in P x and chaos in P y , getting closer and closer to the 
invariant subspaces after each switch and spending longer between switches 
(Figure 8b). The rate of increase of the intervals between switches depends 
on f3, and is about a factor of 1.01 per switch for (5 = 1.10. Decreasing (5 to 
1.088 (Figure 8c, d) results in growth away from cycling chaos and saturation 
to a periodic orbit. 

These calculations required formulating the map (8-11) in terms of the loga- 
rithms of the variables in order to cope with the large (10 1000 ) dynamic range. 
The one place in which accuracy is inevitably lost is in the switch from P x 
to P y (and back), using (10): the C6 1 2/2 1 5s 1 2/o I ^ term is swamped by the order 
one k. As a result of this, the chaotic trajectories start in exactly the same 
way each time the system switches from one invariant subspace to the other. 

In terms of the ODEs, the trajectory entering a neighbourhood of P x close to x- 
rolls shadows the unstable manifold of that equilibrium point (lying inside P x 
and leading to the chaotic set A x ) until the switch to P y . This is in agreement 
with the ODE behaviour shown in Figure 2. 



5 Analysis of the blowout bifurcation 

We briefly review some definitions. If A is a compact flow-invariant subset 
then we define the unstable set 

W U (A) = {x G R 9 : a(x) C A} (17) 
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Fig. 7. Blowout in the map (8-11): compare with Figure 5. (a) Chaotic behaviour 
of x near P x . For different values of (3 (b: /3 = 3.89, c: (3 = 3.87, d: (3 = 3.85), the 
distance from P x as measured by |yo| decays in (b) and grows in (d). The blowout 
bifurcation occurs for parameters near (c). The straight lines show the average 
growth or decay as predicted by (29), averaged over 1000 iterates. 
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(a) (b) 




iterate switch 



Fig. 8. Resonance in the map (8-11): compare with Figures 3 and 6. (a) and (b): 
(3 = 1.100, showing the approach to cycling chaos and the increase in number of 
iterates between switches (p = 1.011). (c) and (d): (3 = 1.088, showing growth away 
from cycling chaos (p = 0.9984) and saturation to a periodic orbit. The dashed lines 
in (c) and (d) show the predicted rate of increase or decrease of the switching times 
from (39). 

and the stable set (or basin of attraction) 

W S (A) = {x E R 9 : u(x) C A} (18) 

where a(x) (resp. oo(x)) is the limit set of a trajectory of the ODE passing 
through x in the limit t — > — oo (resp. oo). We say a compact invariant set A 
is an attractor in the sense of Milnor if 

£{W S {A)) > (19) 

where £(•) is Lebesgue measure on R 9 . It is said to be a minimal Milnor at- 
tractor if there are no proper compact invariant subsets S C A with £(W S (A)\ 
W S (S)) > [15]. 

As shown in [1] , Milnor attractors can occur in a robust manner if the attractor 
lies within an invariant subspace. Suppose P is an invariant subspace and A is 
a compact invariant set in P such that A is a minimal Milnor attractor for the 



18 



flow restricted to P and such that A has a natural ergodic invariant measure 
\i for the flow restricted to P. It is possible to show (under certain additional 
hypotheses [1]) that A as an attractor in the full system if and only if A(fi) < 0, 
where A is the most positive Liapunov exponent for jj, in a direction transverse 
to P (see [5] for a detailed discussion). If we have access to a normal parameter 
[5] such that we can vary the normal dynamics without changing the dynamics 
in P, we can vary A through zero and observe loss of attractiveness of A at 
what Ott and Sommerer have termed a blowout bifurcation [16]. 

5. 1 Set criticality of a blowout bifurcation 

Ott and Sommerer identify two possible scenarios at blowout. At a supercritical 
blowout the attractor bifurcates to a family of attractors displaying on-off 
intermittency [17], with trajectories that come arbitrarily close to A and thus 
linger near A for long times (but with a well-defined mean length of lingering 
or 'laminar phase'). At subcritical blowout there are no nearby attractors after 
loss of stability of A. Note that [3] discusses the question of how to distinguish 
these cases. 

However, what we see in this paper is that there is at least one other possibility 
at blowout bifurcation that is also robust to normal perturbations, namely a 
bifurcation to a cycling state or a robust heteroclinic cycle between chaotic 
invariant sets. This can be seen to be set supercritical but not supercritical, 
in the following sense. 

By reparametrising if necessary we can assume that we have a normal param- 
eter 

c = A(fj.) (20) 

so A undergoes a blowout bifurcation at c = 0. By the argument above, A is 
an attractor only if c < 0. 

If there is a family of minimal Milnor attractors A c , c > such that for all 
e > 



U A C DA (21) 

ce(0,e) 

then we say that the blowout is set supercritical. Note that as discussed in [16] 
and implied by the very word 'blowout', we cannot typically expect the limit 
of the attractors to just be the set A. 
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If in addition A c supports a family of natural ergodic invariant measures fi c 
(c > 0) and — > then we say the blowout is measure supercritical 
or just supercritical (Convergence is the in the C* topology on probability 
measures.) This definition of criticality was used in [3]. 

If the blowout is not set supercritical then we say the blowout is subcritical. 

It seems that one will often get bifurcations in symmetric systems that are 
set supercritical but not supercritical. For example one can generically get a 
blowout from a group orbit of attractors (under a finite group) which yields 
a single attractor that limits onto the whole group orbit. Moreover, in the 
cycling chaos studied here, the attractor after blowout includes not only the 
chaotic invariant sets, but also fixed points involved in the heteroclinic cycle. 
As it is a cycle, it does not possess a natural ergodic invariant measure [21] 
and averages of observables from the system in this state will typically not 
converge but continue to oscillate more and more slowly. 

The blowout scenario described above holds only for variation of normal pa- 
rameters; in general the variation of a parameter will affect both the normal 
dynamics and the dynamics within the invariant subspace. If the dynamics 
in the invariant subspace is chaotic, we can expect to see a large number of 
bifurcations happening within the invariant subspace and these will cause the 
blowout to be spread over an interval of parameter values; there is no reason 
why A(/i x ) should vary continuously with a parameter that varies \i x in a very 
discontinuous manner. 

Nevertheless, for the numerical results presented in Section 3 the dynamics 
in the invariant subspace vary in quite regular way. This is because in our 
system the parameter (3 is a normal parameter for the attractor A x ; for more 
discussion of normal parameters, see [5,8]. 

5.2 Evidence of blowout bifurcation in the ODE model 

By mechanisms described in [20] the pure x and y chaotic dynamics corre- 
sponding to dynamics within P x and P y respectively can become chaotic by 
means of a symmetric global bifurcation that generates Lorenz-like attractors 
approaching the equilibrium solution with full symmetry and x-roll or y-roll 
equilibrium solutions in R x or R y . 

Now suppose there exist chaotic attractors A x and A y contained in P x and P y 
(on average they have more symmetry), From here on we will mostly discuss 
A x but note that as A y is a conjugate attractor, the same will hold for A y . 

These attractors contain a saddle equilibrium e in F (the trivial solution) and 
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so in particular they cannot be Liapunov stable attractors. This is because 
W u (e) H P x must be non-trivial manifold (otherwise e is an attractor; thus 
W u fl P y is also non-trivial as lso{P y ) is conjugate to Iso(P-r). Therefore 

W U (A X ) £ A x (22) 



and so A x cannot be Liapunov stable. However it is possible that A x is an 
attractor in the sense of Milnor; this will imply that the basin will be locally 
riddled in the sense of [5]. 

We assume that A x and A y are minimal Milnor attractors containing e. We 
also assume that they have natural ergodic invariant measures \i x and fi y 
supported on them. We now concentrate on A x . Whether A x is an attractor 
depends on its the spectrum of normal Liapunov exponents. Note that the zero 
Liapunov exponent corresponding to time translation always corresponds to 
a perturbation tangential to the invariant subspace. If 

A(ji x ) < 0, (23) 



where A(/i) is the most positive normal Liapunov exponent for the measure [/,, 
it is possible to show that A x satisfies (19) [5], and hence that A x is a Milnor 
attractor. 

It is comparatively easy to compute A(/i x ) in that case, as the largest trans- 
verse Liapunov exponent of A x in the parameter regime discussed corresponds 
to directions in the Q x direction, where linearised perturbations are described 
by 

y = (v + e(t)-(3x (t) 2 )y . (24) 



Thus we can see that 

AM^pi+(9) flx -P{xl) flx , (25) 



where {f(x))^ = J f(x)dfx(x). From the equation for 9 it is possible to show 
that 



lim 



T T 

f j9dt=\im-Jxtdt (26) 



along bounded trajectories of the ODE. For the given parameter values it is 
possible to approximate (9)^ = —0.03703 and so (xq)^ = 0.03703. Thus 

A(fx x ) = ii - (xl), x - [3(xl), x = 0.12847 - 0.03703/5, (27) 
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implying that the blowout bifurcation occurs at approximately (3 = 3.47, which 
is in good agreement with the simulations (Figure 5b). 

Note that whenever e is hyperbolic and within the attractor A x there exists 
at least one ergodic invariant measure \i e (in particular that one supported 
on e) such that A(/i e ) > 0. In particular, this means that the basins of the 
attractors A x and A y are riddled for all parameter values in our problem! 

Likewise, in the map (8) near P x , the logarithm of y obeys 

log | I log 1 2/o I +logC 2 + <5 2 log \x \, (28) 

where 5 2 is a function of (5. The most positive normal Liapunov exponent in 
this case is then 

A = logC 2 + 5 2 (log|x |), (29) 

where the average is taken over the Lorenz attractor. Averaging over 10 7 
iterates of (7) yields (log \x \) = —5.9233 and, solving (29) for f3, we obtain f3 = 
3.869 at the blowout bifurcation, in agreement with the data in Figure 7(b). 

As A(/j, x ) passes through 0, A x loses stability and becomes a chaotic saddle, and 
in doing so it creates a continuum of connections from A x to the fixed point (y- 
rolls) in R y . These connections are robust to G-equivariant perturbations, as 
2/-rolls are sinks within Q x , and so there is a robust cycle alternating between 
the equilibrium points in R x and R y and the chaotic sets in P x and P y . 

We observed in sections 3 and 4 that this cycling chaos is attractive once it is 
created; we turn to the stability of cycling chaos in the next section. 



6 Analysis of the resonance of cycling chaos 

It is clear from Figures 3, 6 and 8 that the key to understanding the loss 
of stability of cycling chaos lies in obtaining the rate at which trajectories 
approach that state. It is possible to estimate this rate for the map, and we 
use information gained in this calculation to carry out the same estimate in 
the ODEs, and thus are able to obtain the values of f3 at which the loss of 
stability occurs, in the map and in the ODEs. What is remarkable is that, 
once a single average over A x has been computed numerically, the value of (5 
at the bifurcation point can be obtained analytically. 
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6. 1 Rate of approach to cycling chaos 



We suppose that (in the map) the system arrives near P x with given values 
of (xo, 1, yo, 1/2), iterates using (8) n times (n ^> 1), ending up in a state 
(x' , 1, y' , y'2) with \y' Q \ > \x' \. There follows one iterate of (10), leaving the 
system near P y in a state (x'q, x'2, y ', 1) (we are ignoring changes of sign of 
the variables). We need to establish an estimate of (x'q,x'2), a measure of the 
distance from P y , given that (1/0,1/2) were small when the system started close 
to P x . 

Properly, the value of y' after n iterates will depend on the values of x over 
those n iterates, but if n is large, we approximate the detailed history of x 
by its average and obtain 

l°g bo I = log|j/o| +nA e , (30) 
log I2/2I = log I2/2I +nA c , (31) 

where 



Ae = logC 2 + 5 2 (log|x |), (32) 
A c = logC 3 + 5 3 (logk |) (33) 

are Liapunov exponents in the expanding and contracting directions around P x . 
A e is precisely the Liapunov exponent that went through zero at the blowout 
bifurcation in (29). Note that we have effectively approximated the chaotic 
set A x by an equilibrium point. 

The trajectory escapes from the neighbourhood of P x once \y' \ > \x' \; since x 
is typically of order one (compared to the tiny initial value of y ), we assume 
that the escape takes place when \y' \ = 1, so obtaining 

n = -Y l °s\Vol ( 34 ) 
log I2/2I = log I2/2I - ^loglyol, (35) 

with x' and y' both of order one. 
One iterate of (10) now yields 

log I x'o I = log |4I +logC4 + 5 2 log|t/ol +^log|y 2 |, (36) 
log la^'j = log C 5 + 5 3 log I2/0I +<Uog|y 2 |. (37) 



23 



These expressions will be dominated by log 1 2/2 1 once the trajectory is very 
close to the heteroclinic cycle, so, neglecting terms of order one, we obtain 



log I ^ I 

log 141 



x As: 



/log I2/0I 

\ log J/2 



+ 0(1) 



(38) 



for large negative values of log |y | and log 1 3/2 1 - A conjugate map will describe 
the return from P y to P x . One eigenvalue of the matrix is zero because of the 
way we approximated the behaviour near the chaotic set; the other eigenvalue 
is 

P = 5 5 -5 4 ^, (39) 



which we refer to as the stability index. 



The zero eigenvalue forces log \xq\ = (S4/85) log ja^'l, s ° after one iterate of the 
composite map (38), the dynamics will obey 

logKh_ /logjyoh (40) 



log 14V V log 1 2/2 



Clearly if p > 1, log \x \ and log \y \ will tend to —00 and the trajectory will 
asymptote to attracting cycling chaos. Conversely, if p < 1, cycling chaos is 
unstable and trajectories will leave the domain of validity of the approxima- 
tions we have made. We can also deduce from (34) that the number of iterates 
between each switch will increase by a factor of p per switch. 

At the point at which cycling chaos is created (as A e increases through zero), 
we see that p is much greater than one, provided that A c is negative and 64 is 
positive, both of which are true in the examples we have discussed. We deduce 
that cycling chaos is attracting near to its creation at the blowout bifurcation. 



6. 2 Loss of stability of cycling chaos 



From the condition that p = 1 at the loss of stability of the chaotic cycle, we 
determine that this bifurcation occurs in the map at j3 = 1.0896, in agreement 
with the data in Figure 8. 

Returning to the ODEs, we observe that in (39) 64 and 65 are ratios of eigen- 
values of x-rolls (proportional to decay rates of y$ and 2/2), while A c and A c 
are the growth rate of y and the decay rate of y 2 near the chaotic set A x . In 
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the ODEs, the linearisation of (1) about A x yields for the decay rate of 1/2, 
while the growth rate of y is given by (29). Hence we have the stability index 

9 = 65 ~ 64 % (41) 



for the ODEs, where A e is given by (29). Note that #4 and A e are both functions 
of (3. The condition p = 1 is readily solved for (3, and has solution j3 = 1.63. At 
(3 = 1.63 the ODEs are still approaching cycling chaos, with switching times 
increasing by a factor of above 1.02 per switch (see Figure 6). However, the 
ODEs have not yet reached their asymptotic rate of slowing down, resolving 
this small discrepancy. 

On decreasing (3 below 1.63, the stability index p increases above unity and 
the cycling chaos will no longer be attracting. This loss of stability can broadly 
be classed as a resonance of Liapunov exponents. Observe that the resonance 
will be located at different (3 for different invariant measures supported on A XjV 
and so we expect the presence of riddling and associated phenomena found in 
[2] at a resonance of a simpler model displaying cycling chaos. 

We observe that for f3 below the resonance bifurcation, the system exhibits 
behaviour that is numerically indistinguishable from periodic: there appear 
to be a large number of coexisting apparently periodic orbits. We hypothesize 
that the resonance creates a branch of 'approximately periodic attractors', i.e., 
attractors that have a well-defined finite mean period of passage around the 
cycle, going to infinity at the resonance [2]. These might lock onto long periodic 
orbits for progressively smaller j3, as found in the numerical simulations. For 
this example, the approximately periodic attractor branches set supercritically 
from the cycling chaos; however one presumes that subcritical branching is also 
possible. Research is presently in progress on understanding the more detailed 
branching behaviour at this bifurcation. 



7 Discussion 



This study has shown that one possible, apparently generic scenario for loss 
of stability of a chaotic attractor in an invariant subspace on varying a normal 
parameter is as follows: there is a blowout bifurcation that creates an attract- 
ing, robust heteroclinic cycle between chaotic invariant sets (cycling chaos). 
The bifurcation is set supercritical but not supercritical, i.e., the bifurcated 
attractors contain the attractor for the system in the invariant subspace, but 
unlike in a supercritical bifurcation (to an on-off intermittent state) the length 
of laminar phases increases unboundedly along a single trajectory even at a 
finite distance from the blowout. 
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This cycling chaotic state can be modelled well by the network shown in 
Figure 4 although in reality the network is complicated by the facts that (a) 
there are other fixed points contained in the closure of unstable manifolds and 
(b) the fixed points in F x and F y are actually contained in the chaotic sets A x 
and A y rather than being isolated. We suspect this may have the consequence 
that there is no Poincare section to the flow and so the cycle is 'dirty' in 
the terminology of [2]. Nevertheless, the normal Liapunov spectrum of the 
invariant chaotic set seems to determine the attraction or not of the cycle. 

The attracting cycling chaos is observed to lose stability via a mechanism 
that resembles a resonance of eigenvalues in an orbit heteroclinic to equilibria. 
Such a resonance has been seen to occur in special classes of systems with 
skew product structure [2], in analogy to the branching of periodic orbits at 
a homoclinic resonance investigated by [7]. 

Throughout this investigation, it has been necessary to monitor carefully sev- 
eral numerical effects. In particular, for trajectories that display asymptotic 
slowing down characteristic of cycling chaos there will be a point at which 
rounding errors cause the dynamics either to transfer to the invariant sub- 
space, or keep the dynamics a finite distance from the invariant subspace. 
In the context of physical systems there will always be imperfections in the 
system and noise that will destroy the invariant subspaces. Nevertheless the 
perfect symmetry model will be very useful in describing what one expects to 
see in such imperfect situations. 

It still remains to prove rigorously that the observed scenario is generic and 
so of interest to other, less specific models and in particular PDE models of 
which this is a truncation. We could like to emphasise that the behaviour 
we see occurs for a reasonably large region of physically relevant parameters 
in the ODE model and moreover we are unaware of any other classification 
which explains and predicts the observed dynamics to the degree that we have 
done here. In principal cycling chaos can be seen in ODE models down to 
dimension 4, though not smaller than this; thus this dynamics should be seen 
as something that will not be created at a generic bifurcation from a trivial 
state, but rather in a more complicated dynamical regime far from primary 
bifurcation. 

We have discussed a possible route to cycling chaos through a blowout bifur- 
cation, and how cycling chaos might cease to be attracting, both in general 
terms and in the context of a specific model. Our general results ought to be 
applicable to a variety of other problems. In particular, behaviour that might 
be understood in terms of cycling chaos has been seen by Knobloch et al. [12] 
in an ODE model of the dipole-quadrupole interaction in the solar dynamo. 
In their model, there is a weakly broken symmetry between the dipole and 
quadrupole subspaces, and the system switches between the two subspaces, 



26 



favouring one over the other since they are not equivalent. 



Finally, one comment that deserves to be made is that the choice of (3 as 
the parameter allows an important simplification because of this parameter is 
normal for dynamics within P x and P y . One assumes that similar behaviour 
will be observed for non-normal parameters with the exception that the chaos 
in the invariant subspace will be fragile [6] and destroyed by many arbitrarily 
small perturbations; see [8]. 



Acknowledgements 

We acknowledge very interesting conversations with Mike Field, Marty Gol- 
ubitsky and Edgar Knobloch concerning this work. The research of PA was 
partially supported by a Nuffield 'Newly Appointed Science Lecturer' award 
and EPSRC grant GR/K77365. AMR is grateful for support from the Royal 
Astronomical Society. 

References 

[1] J.C. Alexander, I. Kan, J. A. Yorke and Zhiping You. Riddled Basins. Int. 
Journal of Bifurcations and Chaos 2 (1992) 795-813. 

[2] P. Ashwin. Cycles homoclinic to chaotic sets; robustness and resonance. Chaos 
7 (1997) 207-220. 

[3] P. Ashwin, P. Aston and M. Nicol. On the unfolding of a blowout bifurcation. 
Physica D to appear (1997). 

[4] P. Ashwin, J. Buescu and I.N. Stewart. Bubbling of attractors and 
synchronisation of oscillators. Phys. Lett. A 193 (1994) 126-139. 

[5] P. Ashwin, J. Buescu and I.N. Stewart. From attractor to chaotic saddle: a tale 
of transverse instability. Nonlinearity 9 (1996) 703-737. 

[6] E. Barreto, B. Hunt, C. Grebogi and J. Yorke. From high dimensional chaos 
to stable periodic orbits: the structure of parameter space. Phys. Rev. Lett, to 
appear (1997). 

[7] S.-N. Chow, B. Deng and B. Fiedler. Homoclinic bifurcation at resonant 
eigenvalues. J. Dyn. Biff. Eqns. 2 (1990) 177-244. 

[8] E. Covas, P. Ashwin and R. Tavakol. Non-normal parameter blowout bifurcation 
in a truncated dynamo model Preprint, University of Surrey (1997). 

[9] M. Dellnitz, M. Field, M. Golubitsky, A. Hohmann and J. Ma. Cycling chaos. 
IEEE Trans. Circuits and Systems-1 42 (1995) 821-823. 



27 



[10] M. Field. Lectures on bifurcations, dynamics and symmetry. Pitman Research 
Notes in Mathematics, 356, 1996, Pitman. 

[11] J. Guckenheimer and P. Holmes. Structurally stable heteroclinic cycles. Math. 
Proc. Camb. Phil. Soc. 103 (1988) 189-192. 

[12] E. Knobloch, S.M. Tobias and N.O. Weiss. Modulation and symmetry changes 
in stellar dynamos. Mon. Not. Roy. Astr. Soc. to be submitted (1997). 

[13] M. Krupa and I. Melbourne. Asymptotic stability of heteroclinic cycles in 
systems with symmetry Ergod. Th. & Dynam. Sys. 15 (1995) 121-147. 

[14] P.C. Matthews, A.M. Rucklidge, N.O. Weiss and M.R.E. Proctor. The three- 
dimensional development of the shearing instability of convection. Phys. Fluids 
8 (1996) 1350-1352. 

[15] J. Milnor. On the concept of attractor. Commun. Math. Phys. 99 (1985) 
177-195; Comments Commun. Math. Phys. 102 (1985) 517-519. 

[16] E. Ott and J.C. Sommerer. Blowout bifurcations: the occurrence of riddled 
basins and on-off intermittency. Phys. Lett. A 188 (1994) 39-47. 

[17] N. Piatt, E.A. Spiegel and C. Tresser. On-off intermittency; a mechanism for 
bursting. Phys. Rev. Lett. 70 (1993) 279-282. 

[18] W.H. Press, B.P. Flannery, S.A. Teukolsky and W.T. Vetterling Numerical 
Recipes - the Art of Scientific Computing. (Cambridge University Press, 
Cambridge, 1986). 

[19] A.M. Rucklidge and P.C. Matthews. The shearing instability in 
magnetoconvection, in Double-Diffusive Convection, eds. A. Brandt and 
H.J.S. Fernando) (American Geophysical Union, Washington, 1995) pp. 171- 
184. 

[20] A.M. Rucklidge and P.C. Matthews. Analysis of the shearing instability in 
nonlinear convection and magnetoconvection Nonlinearity 9 (1996) 311-351. 

[21] K. Sigmund. Time averages for unpredictable orbits of deterministic systems. 
Annals of Operations Research 37 (1992) 217-228. 



28 



